Equilibrium and Kinetics: Water Confined in Carbon Nanotube as ID Lattice Gas 
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A simple ID lattice gas model is presented, which very well describes the equilibrium and kinetic 
behaviors of water confined in a thin carbon nanotube found in an atomistic molecular dynam- 
ics(MD) simulation [ Nature 414, 188 (2001) ]. The model parameters are corresponding to various 
physical interactions and can be calculated or estimated in statistic mechanics. The roles of every 
interaction in the water filling, emptying and transporting processes are clearly understood. Our 
results indicate that the physical picture of the single-file kinetics is very simple. 



PACS numbers: 47.60. +i, 85.35.Kt, 87.16.Uv 

Transport of molecules through pores on nanometer 
scale is important for biology. Carbon nanotubes have 
been widely studied as pores filled with gases [1,2] and 
liquids [3,4]. Recently, using molecular dynamics(MD) 
simulation, Hummer and coworkers [5] studied the filling 
and transporting of water molecules in a very thin carbon 
nanotube which is a (6, 6) tube with two open ends, 4 A 
radius and 5.5 unit cells (about 13.5 A in length). The in- 
teraction between carbon atoms and water molecules was 
expressed by Lennard-Jones(LJ) potential. They found 
that the LJ parameters (e and a) can sensitively affect 
the filling of water molecules. In their work, two group 
LJ parameters were used: (1) unmodified parameters: 
e = 0.114 kcal/mol, a = 3.2752 A; (2) modified parame- 
ters: e = 0.065 kcal/mol, a = 3.4138 A. For the former, 
it was found that the water molecules completely filled 
the nanotube (the average occupancy number n = 5, the 
corresponding density is about 5 times of that of bulk wa- 
ter, where the radius of the nanotube pore is about 0.8 
A). For the latter, water molecules fluctuated between 
the empty (n = 0) and filled (n = 5) states, and it was 
a two-state like behavior, the lifetimes of the two states 
follow exponential distributions. The high density fill- 
ing seems be surprising in the unmodified system due to 
the hydrophobic property of carbon nanotube, though 
the sharper distribution of the water molecule binding 
energy inside the nanotube was found in the MD simula- 
tion [5] . The sensitive effect of L J potential in the filling 
of water, the two-state like kinetics, and other remark- 
able MD results need to more clearly be understood in 
physics. 

Here, based on a strict statistic mechanics treatment, 
we present a one-dimensional(lD) lattice gas(LG) model 
of water molecules confined in a thin tube with a fi- 
nite length and two open ends. All parameters of this 
model have explicit physical sources, can be calculated 
from interactions. They are: (1) the local excess chemi- 
cal potential Afj, of water in the tube; (2) the free energy 
contribution £ c from LJ potential of carbon atoms; (3) 
the free energy Eh from hydrogen bonds between water 
molecules inside the tube (the rotational free degree of 



the molecules can be included); and (4) the additive free 
energy £ a of the water molecules located at ends of the 
tube, which come from the interaction with the outside 
water molecules. In this model, LJ potential from the 
nanotube wall only affect the parameter £ c . The fitted 
parameters obtained from MD results [5] are in agree- 
ment with the directly calculated or estimated values 
from the corresponding interactions, that is, the obtained 
equilibrium properties are well consistent with the results 
of MD simulation for both the unmodified and modified 
systems. In the LG model, by analyzing transition states 
in which one or two water molecules partially occupied 
the end(s) of the tube, the whole free energy landscape, 
and the main reaction paths can be obtained, therefore, 
the two-state like kinetics of filling and emptying of water 
can be clearly explained. 

Considering a nanotube filled with water molecules, we 
first partition the tube into N cells along its long axis, 
the length of every cell is corresponding to the length of 
hydrogen bond between water molecules. We form our 
model step by step: (1) we consider one water molecule 
moves into a cell from the bulk water, where the cell 
is only a box without any interaction except limiting the 
water molecule in volume Vb . In this process, the changed 
free energy is Afi = —fx w + fx g , where fj, w is the chemical 
potential of bulk water, fi g =t ln(A 3 /Vo), is the chemical 
potential of ideal gas with density 1/Vb, t is a normalized 
temperature. In this paper, the energy unit is set as the 
room temperature 300 K, if no explicit unit is written. 
A = {2TTh 2 /mk B T) 1 ^ 2 w 0.23 A is the de Broglie thermal 
wavelength of water. Actually, A/i is nothing but — (x^j., 
where fi^. is the local excess chemical potential of wa- 
ter inside the nanotube, about —6.87 kcal/mol w —11.6 
from MD simulation [5]. (2) if there is interaction V e (r) 
between the wall of the cell and water molecule, except 
A//, an additive free energy is 



£c = -tM.w I * ffltp[-V e (r)/t]). 

V JVo 



(1) 



For nanotubes, V e (r) is the Lennard-Jones(LJ) poten- 
tial between the water molecule and all carbon atom of 
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the nanotube. Calculation shows that V e (r) is almost 
independent on the azimuth <f> and axis position z in in- 
finite nanotubes [6], For finite nanotubes, the result is 
also right, except a z dependence near the ends of the 
nanotubes. In (6, 6) nanotube, due to its very small ra- 
dius, for both the unmodified and modified LJ parame- 
ters, V e (r) can be well approximated as a box potential 
well, whose widths (about 0.8 A) are in agreement with 
the interior radius of the tube found in MD [5] , and the 
depths of well arc about —8 ~ —9 and —5 <~ —6, respec- 
tively. It is not very difficult to numerically calculate the 
more exact value of £ c , but it is not necessary; (3) if two 
water molecules locate at two near adjacent cells, a hy- 
drogen bond will be formed, the interaction is denoted 
Vff(r). Under the box potential well approximation of 
V e (r), the free energy can be simply written as, 

S H = -t In-L / dri f dr 2 exp[-V ff (n - r 2 )/t], (2) 

where, due to the limit of the length of hydrogen bond 
and the nanotube wall, the two water molecules will oc- 
cupy two cells Vi and Vi (both volumes are Vb), respec- 
tively. The freedom degree r 2 is limited, only some rota- 
tion can be allowed, we have, 

a\ 3 

£ H = E H -t ln(— ) = E H - n g - m, (3) 
Vo 

where E H is the energy of the hydrogen bond. All other 
contributions, including the rotation of the non-spherical 
water molecules , the rotation of water molecules around 
their hydrogen bond, and the possible angle changes of 
the hydrogen bonds, etc., are writen as the factor a (a > 1 
indicates the existence of these freedom degrees) , the cor- 
responding free energy is noted as /j, . While the number 
(n) of water molecules formed a continuous hydrogen- 
bond chain inside the tube is not very large, we suppose 
that fa is independent on n; (4) if a water molecule locate 
at ends of the nanotube, besides £ c (while water is par- 
tially out of the nanotube, even £ c also increase) , there is 
an additive free energy £ a from the interaction between 
the water molecule and the outside water. In fact, £ a is 
also from hydrogen bonds of water molecules, but due to 
the lower density of bulk water, it is higher than £h- 

Based on these considerations, the water molecules in- 
side the nanotube can be treated as a simply ID LG 
model with N sites, including three parameters. Only £ c 
is dependent on LJ parameters, and the dependence can 
be approximately calculated. Hamiltonian of the model 
is, 

N N-l 

U = ^2,£i Ci + £ H ^2 c * c m, ( 4 ) 

i=l i=l 

where £i = A/j + £ c + (S^i + (^jv) £a, <5j,j is Kronecker 
delta, a = 1 (occupied), or (unoccupied). The model is 



equivalent to an Ising model, if we transform Sj = 2cj — 1. 
For a ID Ising model, if the external field H = 0, the 
states with M and — M are symmetry, and the energy 
barrier is very low, where M = ^y S ' . In the LG model, 
the corresponding external field H = — (£, + £h)/2- If 
£ c is very low, H > 0, water molecules can fill the nan- 
otube, while £ c increases (the well depth decreases), H 
decreases. If H <~ 0, there is an equilibrium between 
the filling and emptying states. However, some details 
must be considered as we try to explain the MD results 
of Ref. [5] by using the LG model: (1) the finite sites 
N; (2) the additive energy £ a at ends. Since N is not 
very large, we can directly enumerate all 2^ states. For 
general case, Hodak and Girifalco [7] have written the 
partition function for the LG model with any finite sites 
and the additive energy. 

For unmodified and modified LJ parameters, Hummer 
et al. have given the probability of finding exactly N wa- 
ter molecules inside the nanotube (Figure 2 of Ref. [5]). 
In their MD simulation, N = 5. We can obtain the model 
parameters (£h ~ —6.27, £ a ~ —2.93, and £ c w —4.96) 
by fitting the MD probability data of the modified sys- 
tem. The regenerated free energy profile is very well in 
agreement with the MD probability distribution in the 
modified system as shown in Fig.l. Here, the fitting £ c 
is well consistent with the estimated value from LJ po- 
tential by Eq.(l). The hydrogen bond of water inside 
the nanotube is about 7 kcal/mol s=s 11.8 (Figure 3(a) 
of Ref. [5]), is stronger than their average value in bulk 
water (5 kcal/mol w 8.4). From Eq.(3), we know that 
the free energy contribution /_tj of all rotational freedom 
degrees is about 0.7. The £ a is only half oi£n, the phys- 
ical reason is the lower density of bulk water. Since both 
£ a and £h are independent of LJ potential, their val- 
ues can be used in the unmodified system. For the large 
statistic errors of MD data in the unmodified system, re- 
placing to fit MD results, we directly use that estimated 
£ c rs —8. The obtained free energy profile is shown in 
inset of Fig.l, it is approximately in agreement with the 
MD probability distribution. For equilibrium states, due 
to the large £r, water is favor to the states with a single 
sequentially occupied sites, for example, P(10100) is only 
about 1/500 of P(11000), where P is the probability of 
states. Similarly, water is favor to the states whose end 
sites are occupied. For example, P(01100) is only about 
1/20 of P(11000) due to £ a . For a good approximation, 
we can calculate the partition function only considering 
the states with a single sequence of hydrogen-bond chain 
and at least occupied an end of the tube. We call it as 
single-sequence end occupied approximation (SSEA). In 
a very current MD simulation [8], for short nanotube(13.5 
A) and long nanotube(27 A), Hummer and coworkers 
carried to a longer time MD simulation, by changing LJ 
parameters by varying a control variable 77, e = eirj 2 , 
and a — eri/rj 1 / 6 , where 77 = 1.0 is corresponding to the 
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unmodified LJ parameters. The model can very well ex- 
plain these results. For example, (1) very rare non-SSEA 
events are found; (2) the approximated linear free energy 
profile G{n) for 1 < n < N (Using the SSEA of the LG 
model, G(n) — nAG+const, where AG = A[i+£ c +£h .); 
(3) AG is proportional to r] (Our calculations show that 
the LJ interaction well depth is proportional to r], there- 
fore AG is proportional to rj.). Actually, for explaining 
the equilibrium property found in their MD, Hummer 
and coworkers have formed a free energy formula (Eq.(3) 
of Ref. [8]), which is equivalent to the SSEA of our model. 

The equilibrium property of water filled into nanotube 
has been very well described in the LG model. Be- 
fore investigating the kinetics, we discuss that the states 
with 6 water molecules inside the nanotube, its proba- 
bility is not very small in MD simulation [5,8]. It may 
be corresponding to a state where the ends of the nan- 
otube are partially occupied by two water molecules as 
shown in Fig. 2. The additive free energy of a molecule 
partially occupied is denoted as £' a , which is different 
from £ a , due to both differences of carbon- water inter- 
action and the interaction with the outside water. Since 
water molecules must pass through these partially oc- 
cupied states as their filling or emptying, the free en- 
ergy may be comparable with the kinetic energy bar- 
rier. For the states which n < 6, the effects of the 
partially occupied sites can be ignored, then we have 
£' a = £ a + [G(6) - G(5) -£ H -Afj,- £ c }/2. We obtain 
that £' a is —0.8 and —1.76 for unmodified and modified 
system, respectively. We write that £' a = (a — l)£ c + 5 W , 
where a is a corrected factor of the carbon-water inter- 
action while water molecule partially move out of the 
nanotube, and S w is the contribution from the outside 
water near the end of nanotube. If we suppose that both 
a and 5 W are independent of LJ parameters, we have 
a rj 0.67 and S w rj —3.4. The obtained factor a is well 
in agreement with the calculated result of LJ potential 
while the water molecule is halfly out of the tube, and 
8 W is comparable to £ a . It maybe indicate a real physical 
picture. 

The kinetics of filling and emptying is determined by 
the real reaction paths. The existence of two minimums 
of the free energy profile G(n) is not a sufficient condition 
for the two-state like behavior of filling and emptying. 
The understanding of the whole free energy landscape is 
necessary. Generally, a rate equation can be written as, 
= E b Kb -> a)P b - E fc Ha -» b)P a , where P a is the 
possibility of state a, k(b — > a) = 7 exp[— P(F aj t, — Ft,)] is 
the rate from state b to a, where 7 is a factor, (3 — 1/fcsT, 
F a ,j, is the free energy barrier or saddle point (transition 
state) separating the state b and state a [9,10], F^ is the 
free energy of the initial state b. Obviously, only two con- 
nected states which they have a common boundary can 
transit each other. For the filling, emptying or transport- 
ing of water, there are 2 N states (here we have supposed 
that the outside water molecules be always in equilib- 



rium). The basic motion is that a water molecule jump 
a step inside nanotube or pass through an end of the 
tube. Since the barrier of broken hydrogen bonds is too 
high, the collective motion of a hydrogen-bond chain of- 
ten occure [11], it means all molecules of the chains move 
a basic step at the same time. For any cases, the axis 
positions z of the all moving molecules are good reaction 
coordinates. In the motion of passing an end of nan- 
otube, the free energy barrier locates nearly at the end 
point of the tube (molecule half occupy the end cell), 
so we approximately have the barrier £t ~ £' a , but for 
the motion of water inside of the tube, the barrier is 
very small (w 0) due to the independence of LJ potential 
on z. For the collective motion, the free energy bar- 
rier is the sum of the barrier of every moving molecule. 
Besides this, the facter 7 is different from that of sin- 
gle molecule motion. Actually, 7 is dependent on the 
thermal velocity vt ~ \JksF '/m cx l/^/n, where n is 
the number of molecules in the collective moving chain. 
While the length of a cell I is very small, 7 = vt/1 [12], 
so 7 w 1.5/v 7 ^ ps _1 . 

Consequently, based on the simple estimation, we can 
draw the whole free energy landscape. Obviously, any 
non-SSEA state can fastly transit to the connected SSEA 
states, but the inverse transition is very slow, and the 
barrier between neighbor SSEA states is lower than that 
between the SSEA states and non-SSEA states. For 
example, according to our estimation, from (10000) to 
(11000), the overcome barrier is about AG+£' a —£ a ~ 1.4, 
but from (10000) to an non-SSEA state (10001), the over- 
come energy is about A/j, + £ c + £' a Rj 4.8. Therefore, 
the reaction paths almost do not pass through the non- 
SSEA states in the filling, emptying and transporting 
processes. It indicates that SSEA describes the kinet- 
ics behavior very well. Based on their MD results that 
the probability of non-SSEA states is very low at equi- 
librium, Hummer and coworkers have formed a kinetic 
model only including the reaction between the neighbor 
SSEA states (Eq.(2) of Ref. [8]), and found that it worked 
very well. However, the reason of all reaction paths in- 
cluding non-SSEA states can be neglected is the high free 
energy barriers, not only the high free energy itself of the 
states. 

Now, we can easily explain the two-state like behav- 
ior of filling and emptying: water molecules collectively 
moving a step into or out from the nanotube, the barri- 
ers of the middle states (1 < n < N) are small (about 
1.4 and 1.1, respectively), and the effect of the pre- 
exponent factor is not obvious, so the lifetime of these 
states is very short, the filling and emptying processes 
have a two-state like behavior. Actually, we can cal- 
culate the rate constant of every reaction and the life- 
time of every middle state. For example, we have that 
fc(10000 -> 00000) Rj 0.5, the rate constant of the in- 
verse process is about 0.012, where unit is ps _1 . The 
order of magnitude of the values is in agreement with 
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the estimated results of MD (about 0.3 and 0.009 [8], 
respectively.). 

As this work was being finished, we noticed a work [13] , 
in which Maibaum and Chandler treated whole system 
consisting of the channel of nanotube and the outside 
bulk liquid as a lattice gas by a MC simulation. It was 
shown that there are similar behavior found in MD sim- 
ulation in generic liquid states. Under some approxi- 
mations, they also analytical analyzed the equilibrium 
property of the system. However, for quantitatively ex- 
plaining the behavior of water confined in nanotube, we 
should consider the characters of water, for example, the 
average density of bulk water is only 1 /5 of that of the 
filled state in the tube, then the interaction (£ a ) between 
a water molecule at an end of the tube and the outside 
molecules is weaker than the interaction between two ad- 
jacent inside molecules. For the equilibrium property of 
system, the whole outside bulk water can be thought as 
an aqueous reservoir and contribute an additive free en- 
ergy £ a to the inside molecules at ends of nanotube, but 
for the kinetics, we should identify these states, in which 
different number outside water molecules contact with 
an inside molecule one by one by hydrogen bonds. If we 
denoted the state in which m outside molecules contact 
with n inside molecules one by one as "(m)_[n]" [14], then 
from the state, after a filling or emptying step, the final 
state is "(to — l)_[n + 1]" or "(to + l)_[n — 1]", respec- 
tively. The complete rate chain will consist of all states 
with any to (from to oo). However, since the longer 
outside molecule chain is more instable, in our interest- 
ing time scale, the longer chain will be thought as rapidly 
equilibrating with the surrounding water. Depending on 
the considered largest to, we have different order approx- 
imations of the complete rate equation. In this paper, we 
used the lowest order (0— order) approximation in which 
we supposed that one outside molecule rapidly contacted 
with and departed from the end inside molecule, the equi- 
librium can rapidly reach in our interested time scale. If 
we consider higher order (1— order) approximation, the 
details of energy barrier will be changed, but the differ- 
ences between the two approximations is not very large, 
so we think that the principal character of the kinetics 
of the system has been captured by the simple ID LG 
model. We wish to give more quantitative agreement on 
rate constants by using higher order approximation in 
next work. 
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FIG. 1. For the system with modified LJ parameters (see 
text), the free energy from LG model and —In p(n) from MD 
simulation are plotted, where p(n) is the possibility of finding 
exactly n water molecules inside the nanotube, the simula- 
tion data come from Ref. [5]. The results of the system with 
unmodified LJ parameters are shown at the inset. 
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FIG. 2. n = 6 water molecules are found inside the nan- 
otube, where two water molecules only partially occupy at 
two ends of the nanotube. 
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